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Abstract 

A technique is presented for solving the inverse dynamics and kinematics of 3 degree of 
freedom spatial flexible manipulator. The proposed method finds the joint torques necessary to 
produce a specified end effector motion. Since the inverse dynamic problem in elastic manipulators 
is closely coupled to the inverse kinematic problem, the solution of the first also renders the 
displacements and rotations at any point of the manipulator, including the joints. Furthermore the 
formulation is complete in the sense that it includes all the nonlinear terms due to the large rotation 
of the links. The Timoshenko beam theory is used to model the elastic characteristics, and the 
resulting equations of motion are discretized using the finite element method. An iterative solution 
scheme is proposed that relies on local linearization of the problem. The solution of each 
linearization is carried out in the frequency domain. 

The performance and capabilities of this technique are tested through simulation analysis. 
Results show the potential use of this method for the smooth motion control of space telerobots. 

1. Introduction 

The assumption that robot arms are rigid is not valid when considering light manipulators 
handling heavy payloads with accuracy. In modelling and controlling these manipulators the 
flexibility of their members must be considered. In particular, residual vibrations must be avoided 
when positioning the end-point of a robotic arm. Different techniques, too numerous in fact to be 
mentioned here, have been proposed to achieve this goal. They differ in the mathematical model 
used to describe the dynamics of the system, and by the sensors and control laws used to achieve 
the overall motion of the arm while controlling the tip vibrations [1-5] .These results show either an 
initial motion of the tip in the opposite direction of the overall motion (undershoot) or an overshoot 
of the set-point and residual vibrations at the end of the trajectory. The initial undershoot is 
characteristic of the step response of a non-minimum phase system, such as a flexible link [1-5]. 
By relying entirely on feedback schemes, the mentioned techniques implicitly considered causal 
controls. 

The authors introduced a direct approach for the solution of the inverse dynamics of a single- 
link [6,7], and multi-link [8] flexible manipulators. The proposed method finds the joint torques 
necessary to produce a desired end effector motion. The new technique for inverse dynamics has 
been validated by feed-forwarding the computed torques in numerical simulations [6-8], and in 
experiments [7,8] for both single and multi-link planar flexible robots, yielding excellent results. 
For a generic acceleration profile the solution of the inverse dynamics yields a torque that must be 
applied before the end-effector starts actually moving. The inverse dynamics provides a unique 
non-causal control law [9] which is capable of tracking a specified trajectory with no deflections 
from it. In particular, it is possible to track any Fourier transformable acceleration profile with zero 
initial undershoot or overshoot, and this therefore contradicts the belief that a non-minimum phase 
system can not be inverted for any given trajectory [10]. This only applies if causal solutions are 
sought, as is the case in Ref. [10]. 

In this paper that formulation is extended to compute the joint torques (and joint angles) that 
position the end-effector of a 3 degree of freedom spatial flexible manipulator along prescribed 
trajectories. The new formulation includes all the nonlinear Coriolis and centrifugal effects. 
Timoshenko beam theory is used to model the elastic behavior of the system and the resulting 



equations of motion are discretized using the finite element method. An iterative solution scheme is 
proposed for finding the desired joint torques, where the solution of each linearization is carried 
out in the frequency domain. A simulation is performed that demonstrates the capabilities of the 
present formulation in handling the problems of vibration control and trajectory tracking of spatial 
manipulators. The technique is therefore of direct use for open-loop control and it also provides the 
basis for closed-loop control algorithms. 


2. Kinematics 

Let us consider the spatial manipulator depicted in Fig. 1 that has three degrees of freedom 
and three revolute joints. The objective is to move the end effector along a given trajectory without 
overshoot and residual elastic oscillations of the tip. These oscillations are due mainly to the 
transverse elastic displacements of the links. We consider the case in which a general three 
dimensional motion is decomposed into three parts as shown in Fig. 2. The first and third parts are 
planar motions and involve the rotations about the Z2 and Z3 axis. The second part of the motion is 
composed of a rotation about the Z\ axis. 

Kinematics of the first and third motion 

The first and third parts of the motion are planar. Consider the individual flexible link 
depicted in Fig. 3 of total length Lthat forms part of the spatial robot. A point P at a distance x 

from the center of the hub has undergone elastic deflections of value u x and Uy , and rotation 9. 

These are defined with respect to a nominal position characterized by the moving frame (ej ef) that 

rotates at a specified nominal angular velocity and acceleration co n and a n , respectively. This 
definition of the motions with respect to the nominal frame permits the linearization of the problem 
from the outset. 

As a consequence of the elastic deflections and rotating nominal motion, the point P is 
subjected to a total linear acceleration a p and angular acceleration a p . The acceleration of the point 
P can be set in terms of the nominal translational and angular accelerations at the hub, a n and a n , 

angular velocity co n at the hub, and the relative velocity v re i and acceleration a re i of the point P. The 
latter are due to the elastic deflections with respect to the moving frame. The kinematics of the 
motion can be established as follows: 

a p = co n x {co n xr) + a n xr + 2 co n x v re i + «„ + a rei (1) 

a P = a n+V 

where r = {x + uj e 2 (t) + u y e 2 (t) 

The components of the relative velocity are w*. and Uy, and those of the relative acceleration 

are \ and Performing the vectorial operations involved in Eq. (1) the following components 
of the accelerations are obtained: 
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Kinematics of the second motion 

Fig. 4 shows the deformation of the manipulator in this second part of the motion which 
consists in a rotation about the global y axis. The position vector for the point P now becomes. 
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( 3 ) 


r — (x + u x ) €j+ (y + Uy ) ^2 + e 3 

Note that the moving frame is now (ej e$) which rotates about the y axis, hence the vectors 
CD n and a n contain only one non- zero component. Again performing the vectorial operations 
involved in Eq. (1) the following six components of the accelerations are obtained. 

a x = -cc& x - co£ u x + a n u x + 2 C0nU z + ti x , a x ~ 0 X 

a y = Uy , a y = a n + 6y (4) 

a z = - cc£ u z - oc n x - a n u x - 2 co n ti x + U z , a z = 0 Z 

Note that this part of the motion has a full three dimensional behavior, and that in this case 
there is coupling between the x and z components of the motion through the elastic deformations. 


3. Equations of Motion 

Once the displacement field has been defined, the principle of virtual displacements can be 
used to directly set up the dynamic equations of motion. We rely on the Timoshenko beam theory, 
which includes the effects of shear deformation and rotatory inertia. For the general three 
dimensional case the equations become: 
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where rj x , T\y and ri z are the radii of gyration of the section, m the mass per unit length, I the 
moment of inertia, A the area, E the Young modulus, G the shear modulus and shear coefficient k. 

A tip masses are represented by M t , and the hub inertias by 1^. 8u x , 8u,y> 6u z , 86 x , 86 y and 86 z 
represent a set of virtual elastic displacements. T=(Tq, Tj, T 2 ) are the unknown torques to be 
applied at the different joints so that the prescribed end-effector motion is obtained. The equations 
of motion for the planar case (motions 1 and 3) can be obtained as a particular case of Eq. (5), and 
they are described in detail in [8]. Note that the accelerations at the joints will have two 

components a n and*(9. The first are the nominal accelerations and the second are the acceleration 
due to the elastic deflections. 

The displacement field of Eq (5) can be discretized using the finite element method. A 
complete set of interpolation functions are defined within each element: 


n 

Ujx.l) = ‘ 

J i=l 


H.(x) u.(t) and 


Bfx,t) = 2 d H{x) d(t) 


( 6 ) 
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where Hi are interpolation functions whose order depends on the number of nodes, n, in the 

elements; u l j and 6 indicate the components j of the nodal deflections. Substituting Eq (4) and (6) 
in the virtual work expression, Eq (5), and following standard procedures for the formation and 
assemblage of element matrices [11], the equations of motion of the manipulator may be expressed 
by a set of time varying differential equations: 


M v+ [ C + C c (cc n ,co n )] v + [K+K c (a n ,ty n )] v = T - F(a n ,co n ) 


M and K are the conventional finite element mass and stiffness matrices, respectively; C c and K c 
are the time varying Coriolis and centrifugal stiffness matrices, respectively. The matrix C has 
been added to represent the internal viscous damping of the material. The vector T contains the 
unknown torques at the joints. F the known forces produced by the rotating frame effect. 

The set of finite element equations (7) may be partitioned as follows: 





where 6 ^ are the elastic rotations at the joints, v t represent the elastic deflections ( x,y,z ) at the 
end-effector, and v t are all the other internal finite element elastic degrees of freedom of the 
manipulator. The force vector F is partitioned in the same manner. 

4. Inverse Dynamics 

The essence of the inverse dynamics is to find the torques T at the robot joints that will cause 
the end of the manipulator to move along a specified nominal trajectory. The problem can be 
quantitatively stated as the finding of the torques T in Eq. (8) under the condition that elastic 
deflections at the end-effector v t (t) be zero. 


Planar motions 

The first and third motions are planar. In this case the best way to solve the inverse dynamics 
is by establishing a recursive procedure as explained in detail in reference [8]. This procedure 
consists in finding first the torque corresponding to the last link under the condition that the elastic 
displacement at the tip is zero, the next step is to compute the reactions at the hub which will be 
transmitted to the next link in the chain. These reactions may be obtained simply by equilibrium 
considerations. The procedure continues with the next link in the chain in the same manner as 
before with the reaction forces included. Note that this process is conceptually similar to the 
recursive Newton-Euler scheme for inverse dynamics of rigid robots. The method continues with 
the rest of the links until the torque on the first link is determined. This way of proceeding assures 
that the end of each link moves along the desired nominal trajectory without oscillating from it. 
Once the torques T have been obtained, the motions at any point of the links specified by v,- or the 

angles 0 h (inverse kinematics) follow from Eq (8) by a direct analysis. 


Spatial Motion 

In the case of three dimensional manipulators with elastic properties in all directions, there 
are elastic displacements that can not be controlled by the corresponding joint torque. These elastic 
displacements will influence the motions of subsequent links, introducing perturbations in their 
nominal motion. The nominal position of each link will be modified by the elastic displacement at 
the end of the previous link. The recursive procedure mentioned above now requires an iteration 
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process that will account for the displacement corrections starting with the last link. For the spatial 
motion at hand the manipulator does not change its configuration (except for the global rotation), 
and rather than using the recursive procedure, it is easier to find the solution for the three joint 
torques T directly from Eq. (8) under the condition that three elastic deflections at the end-effector 
v t ( t) be zero. 

The forcing terms in the right hand side of Eq. (6) depend on the velocities and acceleration 
profiles of the nominal motion of the link. These will have to be Fourier transformable for the 
problem to have a solution. Again once the torques T have been found, the motions at any point of 

the robot specified by v t - and the joint angles 6 h (inverse kinematics) can be readily obtained from 
Eq.(8) by a direct analysis. It becomes obvious at this point how the inverse dynamics and 
kinematics are closely coupled in the analysis of flexible robots.The torques T can be found in the 
frequency domain by means of an iterative numerical procedure similar to that used by Bayo et al 
[7-8]. In order to set up the iterative process Eq. (8) may be restated as follows: 

Mv+Cv + Kv = IT - F - C c (« n ,ct) n )v - K c (a n , co n )\ (9) 

where all the time invariant terms have been left in the L.H.S. of the equation and the time variant 
ones have been moved to the R.H.S. The vector l contain unit values for the degrees of freedom 
corresponding to the joint rotations and zero for the rest, as shown in Eq. (8). 

The iteration process is initiated by solving Eq.(9) for the unknown torques T in the 
frequency domain (see [7-8] for a detailed description of this process). The first iteration is done in 
the absence of the last two terms involving C c and Kc in the RHS. The first iteration will yield a 

displacement vector v x (t) which in turn will be used to compute K c and C c , and the last two terms 
in Eq. (9). The process is then repeated with the new force vector under the constraint that 
v t (t)=0. The iteration procedure may be stopped when the norm 1171 - TV ill for the solution of two 
consecutive iterations, is smaller than a prescribed tolerance. 

5. Simulation Analysis 

Let us consider as an example a manipulator with the following characteristics: 

First link: L=0.2 m, A=1.90xl0- 4 m 2 , Ii=I 2 = 3.5x10-9 m 4 , 

J=2.0xl0-8m 4 , M t =2.0 Kg, I h =0.002 m 4 

Second link: L=0.6 m, A=1.20xl0- 4 m 2 , Ii=I 2 = 2.2864xlO- 10 m 4 , 

J=4.57xlO- 10 m 4 , M t =1.0 Kg, I h =0.0012 m 4 

Third link: L=0.6 m, A=0.40xl0~ 4 m 2 , I]=I 2 = 8.4683xl0 -12 m 4 , 

J=1.69xl0" 11 m 4 , M t =0.15 Kg, I h =0.00048 m 4 

They share the following properties: E=7.11xl0 10 N/m 2 , mass density = 2715 Kg/m 3 , shear 
coefficient k=l/2 and a damping ratio of 0.002. The first frequencies of the bending modes in 

rd/sec of the second and third links under pinned- free conditions are: for the second link co^O 

(rigid body mode), 0)2=157.53, 0)3=487.55; and for the third link 0^=0, 0)2=54.73, 0)3=181.98 

and 0) 4 =389.93. The natural frequencies of first link are so high that it is considered as being rigid. 
The manipulator is modeled with 1 1 finite elements with a total of 67 degrees of freedom. 

The first and third motions consist in straight-line tip trajectories (each one of 1 second of 
duration) generated according the acceleration profile shown in Fig. 5 which corresponds to a 
Gaussian velocity profile [12], The second motion consists in a circular tip trajectory whose 
acceleration is shown in Fig. 5 also and which last 3 seconds. The solution of the inverse 
dynamics is not limited by the type of trajectories chosen, as long as their acceleration profiles are 
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Fourier transformable [9]. As shown in [7] and [12], other type of trajectories can be perfectly 
used. However the trajectory corresponding to a Gaussian velocity profile is selected because it 
diminishes the high frequency content of the input function. 

The torques that need to be applied at each actuator to yield the required tip trajectory are 
calculated according to the procedures described above. As examples Fig. 6 shows the torque T 2 
(flexible torque) for the first motion, compared to the torque needed to produce the same tip motion 
on an infinitely rigid link of the same mass (rigid torque). Fig. 7a and 7b show the torques Tq and 
T 2 for the second motion. It is important to note that in order to track the desired trajectory the 
flexible torques need to be applied before the tip actually moves, as also noted in [7-8]. In 
addition, it may be seen also that the flexible torques have a smaller peak value, and are more 
spread over time than the rigid one. As a result of the inverse kinematics, Fig. 8 compares the 
nominal angle at joint 2 with the calculated one for the first part of the motion. These could 
obviously be used in a collocated feed-back control scheme. 

All the calculations were performed using a SUN 160 workstation. The total time required 
for the computation of these torques for all the three motions was approximately 20 seconds. The 
computational burden is mostly due to the solution of the complex sets of equations and the use of 
the fast Fourier transform. The computed torques produce the desired tip trajectory without 
overshoot or residual oscillations. A forward dynamic analysis and graphics simulation are 
performed yielding the results shown in Fig. 9 and 10. It may be observed how both links bend 
along their motion. However, their ends follow the specified nominal path. Finally Fig. 12 
illustrates the motion of the end-effector when the rigid torque at joint 0, corresponding to the 
second motion, is input to the system. Although the rigid torque is smooth the tip still vibrates, and 
this points out the fact that smooth trajectories or torques are not enough to eliminate residual 
vibrations. 

6. Conclusions and Acknowledgement 

A procedure has been presented for the inverse dynamics and kinematics of a three degree of 
freedom spatial manipulator. The proposed algorithm calculates the torques required to move the 
end effector through a specified trajectory while avoiding tip oscillations. The tip trajectory 
acceleration profile should be Fourier transformable. The technique is demonstrated through 
simulation yielding very good results. The algorithm converges very rapidly in all the cases 
considered. 

The use of the fast Fourier transform, solution of complex equations and the necessity to 
iterate to obtain the required solution preclude this method from being used in real time control. 
However, it can still be used in open-loop control and other close-loop control strategies. 
Currently, under investigation are time domain algorithms to be implemented in multi-processors 
computer architectures which will substantially reduce the computational time. 

The authors would like to acknowledge the support of this work by NATO, grant number 
0877-87, and by the National Science Foundation, grant number 8421415, through the Center for 
Robotics Systems in Microelectronics. 
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5. Acceleration profiles for the three motions. 



6. Joint 2. First part of the motion. 

Rigid torque (dotted) vs. flexible torque (solid). 
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7a Joint 0. Second part of the motion. 

Rigid torque (dotted) vs. flexible torque (solid). 
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7b Joint 2. Second part of the motion. 

Rigid torque (dotted) vs. flexible torque (solid). 
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8 Joint 2. First part of the motion. 

Nominal (dotted) vs. calculated (solid) joint angles. 



9. Simulation of the first part of the motion. 
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10. Simulation of the second part of the motion. 


1 1 . Vibration at the end-effector when using the rigid 
torque. Second part of the motion. 






